Sensors 2014, 14, 8283-8304; doi:10.3390/sl40508283 



OPEN ACCESS 



sensors 

ISSN 1424-8220 

www.mdpi.com/journal/sensors 

Article 

Dynamic Compressed HRRP Generation for Random 
Stepped-Frequency Radar Based on Complex- Valued 
Fast Sequential Homotopy 

Peng You, Zhen Liu Hongqiang Wang, Xizhang Wei and Xiang Li 

School of Electronic Science and Engineering, National University of Defense Technology, 
Changsha 410073, China; E-Mails: ypnudt@126.com (P.Y.); oliverwhq@vip.tom.com (H.W.); 
liwerer@nudt.edu.cn (X.W.); lixiang01@vip.sina.com (X.L.) 

* Author to whom correspondence should be addressed; E-Mail: zhen_liu@nudt.edu.cn; 
Tel/Fax: +86-731-8457-3493. 

Received: 8 February 2014; in revised form: 24 April 2014 /Accepted: 5 May 2014 / 
Published: 8 May 2014 



Abstract: Compressed sensing has been applied to achieve high resolution range profiles 
(HRRPs) using a stepped-frequency radar. In this new scheme, much fewer pulses are 
required to recover the target's strong scattering centers, which can greatly reduce the 
coherent processing interval (CPI) and improve the anti-jamming capability. For practical 
applications, however, the required number of pulses is difficult to determine in advance and 
any reduction of the transmitted pulses is attractive. In this paper, a dynamic compressed 
sensing strategy for HRRP generation is proposed, in which the estimated HRRP is updated 
with sequentially transmitted and received pulses until the proper stopping rules are satisfied. To 
efficiently implement the sequential update, a complex-valued fast sequential homotopy 
(CV-FSH) algorithm is developed based on group sparse recovery. This algorithm performs 
as an efficient recursive procedure of sparse recovery, thus avoiding solving a new 
optimization problem from scratch. Furthermore, the proper stopping rules are presented 
according to the special characteristics of HRRP. Therefore, the optimal number of pulses 
required in each CPI can be sought adapting to the echo signal. The results using simulated 
and real data show the effectiveness of the proposed approach and demonstrate that the 
established dynamic strategy is more suitable for uncooperative targets. 

Keywords: high resolution range profile (HRRP); random stepped-frequency; dynamic 
compressed sensing; fast sequential homotopy; recursive sparse recovery; group sparse 
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Acronyms: 



HRRP 


High range resolution profile 


CV-FSH 


Complex-valued fast sequential homotopy 


RV-FSH 


Real- valued FSH 


CPI 


Coherent processing interval 


PRF 


Pulse repetition frequency 


PRI 


Pulse repetition interval 


LSFR 


Linear stepped-frequency radar 


SSFR 


Sparse stepped-frequency radar 


RSFR 


Random stepped-frequency radar 


CS 


Compressed sensing 


RIP 


Restricted isometry property 


SCS 


Sequential compressed sensing 


SNR 


Signal-to-noise ratio 


T/R 


Transmission and reception 


BPDN 


Basis pursuit denoising 


LOS 


Line of sight 


CVX 


Convex 



1. Introduction 

Radar can be used to achieve high resolution images of targets. High range resolution profiles 
(HRRPs) are easier to obtain compared to their high dimensional counterparts and can be used in many 
applications including target classification and recognition, etc. Linear stepped-frequency radar (LSFR) 
is an effective way to synthesize high resolution range profiles of targets, which can be implemented by 
a frequency-agile transmitter with a narrowband receiver [1]. Thus, LSFR can be applied to complement 
the imaging capability for narrowband radar systems. In order to synthesize a large bandwidth for LSFR 
with narrowband pulses, one must make an appropriate tradeoff between pulse repetition frequency 
(PRF) and coherent processing interval (CPI) from a traditional view. On the one hand, in the cases of 
fixed low and median PRFs, a large number of pulses should be transmitted to cover the whole 
frequency band, thus resulting in long observation times. However, modern radars often have other 
missions besides imaging, such as detecting and tracking, which will limit the observation time for each 
target. Furthermore, for long CPIs, both radial- and micro-motions of the target will seriously distort the 
pulse coherence and thus degrade the performance of HRRP generation with LSFR. On the other hand, 
when the CPI is fixed, high PRF is usually required for LSFR to achieve enough coherent pulses, which 
will induce range ambiguity and blind zones for farther targets. 

One effective way to achieve large bandwidth with low PRF during short CPI is to employ the sparse 
stepped-frequency radar (SSFR) [2], in which there are far fewer pulses in the train due to the frequency 
band vacancy. For SSFR, there are mainly two sparse patterns, Le,, the periodic pattern and the random 
pattern, both of which will degrade the HRRP quality if we directly apply traditional pulse compression 
techniques. For the periodic pattern combined with Fourier transform, grating lobes will appear to cause 
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aliasing for HRRPs, and for the random pattern, the correlation processing suffers from high sidelobe 
pedestals, which will decrease the peak to sidelobe ratio. 

In order to reduce the grating lobes and sidelobe pedestal caused by missing frequency bands in SSFR, 
high quality HRRPs are achieved in [2,3] by applying the sequence CLEAN method, which is 
essentially a simple sparse-seeking approach and can be viewed as an implementation of matching 
pursuit [4]. More recently, owing to the advances in optimal recovery from sparsity constraints [5], 
especially the emergence of compressed sensing (CS) theory [6,7], exploiting sparsity in high resolution 
processing for SSFR has been attracting growing attention. For example, the sparsity-driven HRRP 
synthesis based on Bayesian CS (BCS) is proposed in [8] for the SSFR with periodic patterns. In [9,10] 
the CS theory is also applied to form HRRPs of moving targets using the SSFR with random pattern. 
However, these are only pilot studies and focus on finding the sparse HRRP for a given set of echoes. 
For practical applications, the required number of pulses is difficult to determine in advance under 
various scenarios and the reduction of the transmitted pulses is attractive. A meaningful solution is that 
we transmit the HRRP pulses and update the estimated HRRP with sequentially received pulses until the 
proper stopping rules are satisfied. In this case, the existing static recovery procedures always suffer 
from low efficiency due to having to solve a new inverse problem from scratch when the echo pulses are 
available sequentially. 

In this paper, we propose a novel algorithm to quickly update the estimated HRRP without solving a 
new optimization problem from scratch as well as establishing the proper stopping rules. Due to the 
linearity of frequency, LSFR as well as SSFR with periodic pattern have the "diagonal ridge" ambiguity 
function, which suffers from serious delay-Doppler coupling [11]. Therefore, herein we only focus on 
SSFR with random patterns, which can also be called random stepped-frequency radar (RSFR) [12]. 
Besides delay-Doppler decoupling, excellent resistance to range ambiguity and electronic countermeasures 
can also be achieved by the randomly transmitted frequencies in RSFR, which even satisfy the request of 
randomness for the restricted isometry property (RIP) condition [13] and guarantee that HRRP can be 
recovered exactly with high probability [10]. The main contributions of this paper are thus in the 
following three aspects: 

(1) A dynamic compressed sensing strategy for HRRP generation is provided. In this strategy, the 
estimated HRRP is updated with sequentially transmitted and received pulses. Furthermore, the number 
of required pulses within each CPI is determined adaptively for attaining an acceptable HRRP in various 
scenarios, e.g., different targets and/or different target aspects. 

(2) A complex-valued fast sequential homotopy (CV-FSH) algorithm is developed to efficiently 
implement sequential update. This approach is an extension to the fast sequential homotopy (FSH) 
proposed in [14], which has been proved more effective for sequential solution in the noisy settings than 
its counterparts [14]. It should be noted that the extension is not straightforward. We first reformulate the 
complex-valued sparse recovery as a group sparse recovery model. Then the fast sequential update for 
group sparse recovery is completed by resorting to the homotopy technique. 

(3) Proper stopping rules for the dynamic approach are given based on the special characteristics of 
HRRP. With these rules, we can seek the optimal number of pulses required in each CPI adapting to the 
echo signal. In the experiments, simulated and real data are used to test the proposed approach. The 
results show its effectiveness for both stationary and moving targets. 
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The paper is organized as follows: Section 2 describes the signal model of RSFR and the sequential 
compressed sensing (SCS) based dynamic strategy for HRRP synthesis is derived. In Section 3, the CV-FSH 
algorithm is developed to implement fast sequential processing of SCS and some proper stopping rules are 
given for the dynamic approach in Section 4. Section 5 presents some results achieved by using both 
simulated and chamber measured data, which validate the effectiveness of the proposed algorithm. The last 
section gives the conclusion. In this paper, the operations of transposition and conjugate transposition are 
denoted by superscripts T and //, respectively. Operators Re(-) and Im(-) select the real part and imaginary 

part of the argument, respectively. ||x||^ denotes the /^-norm operation of x , and |-| denotes the 
absolute operation. 

2. Dynamic High Range Resolution Profile (HRRP) Synthesis in Random Stepped-Frequency 
Radar (RSFR) 

2.7. Echo Signal Model of RSFR 

As shown in Figure 1, the pulse train signal used in RSFR is composed of M pulses, where M is to 
be determined according to various scenarios. The carrier frequency of the mth pulse isF^ = /c +/m ' 
where is the fundamental carrier frequency and = C^Af is randomly distributed in the bandwidth 
B with a basic frequency step A/ . is an integer randomly selected from {0,1,..., A/^-l} , N = B/Af . 
The pulse repetition interval (PRI) is and the pulse width is T . 

Figure 1. Pulse train signal model of Random stepped-frequency radar (RSFR). 




Assume that the extended rigid target has K scattering centers projected on the radar line of sight 
(LOS) and that the aspect of the target with respect to radar remains unchanged during the coherent 
processing interval (CPI) for HRRP synthesis. The intensity of the kth scattering center is and its 

initial distance apart from the radar is R^j^ . When the target moves with a constant radial velocity of v 

(positive when towards the radar), then, after down-conversion and low-pass filtering, the echo signals 
can be represented as: 



k=l 



t-mT-T 12-t, 



exp(-j2nFj,) 



(1) 
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, m = 0,l,...,M-l , is the echo envelop and 

else 



tj^ = — = t^j^ , c is the velocity of light. Sampling the signal at = mT^ ( = 2^ I c , 

c c 

and /?J represents the sampling distance nearest the target), then the output signal can be expressed as: 

K 

sXm\ = J^Al exp[-j2nFJt,, -2vt^/c)] (2) 

k=l 

Assume that the unambiguous delay (i.e., the PRI) T is divided into N (i.e., the number of range cells, 
which is larger than M) points as T =tQ+nAt with At = l/B and = 0,1,..., and the velocity 
domain of interest Vu is divided into L points as =v^+ lAv , where is the initial target radial velocity 
and / = 0,1,...,L-1 . Then, by considering noisy cases. Equation (2) can be expressed as: 

s,=U,,a+n (3) 

where U^^ is an MxN random matrix, whose elements are 
UiQ[m,n] = cxp(-j2TiFJ'Jcxp(j4TiFjJ^iJc) with Viq=v . The vector a =[aQ,a^,...,a^_^f just 

represents HRRP, in which there are about K non-zero elements, n represents white zero-mean 
measurement noise. 

2.2. SCS-Based HRRP Synthesis for RSFR 

Obviously, Equation (3) can be regarded as a typical underdetermined linear system, and the 
randomness of carrier frequencies ensures that the equivalent observation matrix is a randomly sampled 
sub-Fourier matrix, which could meet the RIP condition [15]. Moreover, if the signal-to-noise ratio 
(SNR) is sufficiently high, and the number of main scattering centers projected on the radar line of sight 
is sufficiently smaller than the number of range cells to guarantee the following sparse condition [13]: 

M>o[K'(logNy) (4) 

then Equation (3) can be regarded as a typical CS model in the complex- valued domain and the HRRP 
can be obtained by solving the following basis pursuit denoising (BPDN) problem: 

• 1 II TJ l|2 ^ II II 

a =argmin — lls*^ +^lrlli (5) 

^ 2 

where }i is a regularization parameter. 

In practice, however, the number of main scattering centers is always changing for different targets or 
the same target with various attitudes, thus, the number of pulses is difficult to choose a priori. As the 
observations are available in sequence for RSFR, herein we can assume that the number of pulses in each 
CPI is bounded in a proper scope and then the echo signal can be more accurately modeled as an SCS 
problem [16]. Unlike traditional scheme which is performed after all the echo signals are obtained, here 
we can generate the HRRPs sequentially with echo pulses from the minimum number to the maximum 
number. In this situation, the interval time between each pulse can be exploited for fast sequential 
processing. The dynamic algorithm for HRRP synthesis as well as motion compensation is shown in 
Figure 2 and listed as the following steps. 
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Algorithm 1 Dynamic HRRP synthesis and motion compensation 



Step 1: Initialize the number of pulses as the minimum number, i.e., M =M^^ . Assume the 

unambiguous velocity achieved directly from tracking system or by resolution of any Doppler ambiguity 
(depending on the tracking scheme) is v and the velocity measurement accuracy is ±5v , then the 
velocity domain of interest can be set as v^=[v- 5v, v + 5v] , the initial velocity is: 

v^=v-5v (6) 

and the number of velocity cells is L = 25v/ Av , where Av is the velocity step size. 

Step 2: The HRRP in the /th velocity cell is figured out by resolving the following BPDN problem: 

1 II i|2 II II 

«L =argmin-^^^-t/>^ (7) 

where t/^ is an MxA/^ random matrix with elements Ul^[m,n] = cxp(-j27iF^T^)Gxp(j47iFj^Vi/c) 
and X.^ is the corresponding regularization parameter. 

Step 3: The minimum /i-norm criterion is applied for motion compensation by calculating the /i-norm 
of every and searching the minimum to get: 

(a^,/^) = argmin||a^|| .gx 

where is the estimated HRRP for the M pulses and is the estimation of the corresponding 
velocity cell. 

Step 4: When the (M + 1) the sample of echo signal is obtained, we have s^^^^ = 



Set U'^,,= 



u 



M 

I 

M+l 



s^[M+l]_ 

with ii^^i =[m^^J0],...,m^^i[A^-1]] and 

u^M^^[n] = Qxp(-j2nF^^;rjQxp(j4nF^J^^yi/c). Take «^ , and ^ as the "warm-starts" and 
re-estimate a^^^ , a^^^ and /^^^ by resolving the following two problems: 

1 II i|2 II II 

«Li =argmin-|^^^^i -f^M+i<+i|, (9) 



(10) 



Step 5: If the estimated HRRP and velocity cell for the M+l pulses are sufficiently similar to these for 
the M pulses or the number of pulses reaches the maximum number, i.e., M +1 = , then with high 

probability the a^^^ obtained in Step 4 is the HRRP we need and the estimation corresponding to /^^^ is 
just the target velocity v^^^ = K , else enter a new transmission and reception (T/R), and set 

'■M+l 

M=M+1, s^^ =s^^^,, t/^ =t/Ln =^M+n 4 = , then go to Step 4. 



The fast sequential processing algorithm for update is developed in Section 3. In Section 4, some 
stopping rules are discussed and the proper stopping rules are given. 
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Figure 2. Flow chart of the dynamic HRRP generation processing for RSFR. 
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3. Fast Sequential Processing of SCS 

Although the HRRPs of the target can be obtained well with the dynamic algorithm proposed in the 
last section, there are still two important issues from practical consideration. One is how to avoid solving 
a new optimization problem from scratch when the new echo pulses are available sequentially. The other 
is how to design proper stopping rules for attaining acceptable results in various scenarios, which can 
determine the minimum number of pulses required and make the CPI the shortest. Therefore, in this 
section we first assume that the BPDN (7) has already been solved by using the homotopy technique, 
which is computationally efficient and has a simple regularization parameter setting [17], and then the 
FSH approach [14] is extended to resolve (9) more efficiently. In FSH, homotopy algorithm is modified 
to update the solution for BPDN when a new measurement in the real-valued domain is added to the 
system sequentially. In Subsection 3.1, the complex-valued sparse recovery is transformed into a 
real- valued group sparse recovery model. The quick update of group sparse recovery solution based on 
the homotopy technique is presented in Subsection 3.2. 

3.1. Group Sparse Recovery Model for Complex-Valued Sparse Recovery 



In order to make the deduction simpler, we rewrite the abstract signal model in the complex-valued 
domain as well as the BPDN problem as: 

s=Ua+n 



a = argmin — ||?7(a -5* 

« 2 ^ 



(11) 
(12) 



where JJ gC^""^ is the equivalent sensing matrix, s'^,wgC^''^ are vectors of the measurements and 
complex noise, respectively, a g C^""^ is a sparse vector and X. is a regularization parameter. As shown 
in [18], in order to use the sparse recovery techniques in the real-valued domain, the complex- valued 
model should be firstly transferred to the real- valued one by defining: 

"Re(t/),-Im(l/) 
Im(l/),Re(t/) 



"Re(s,)" 




"Re(a)" 




"Re(n)" 


_Im(5,) 


,x = 


Im(a) 




Im(n) 



(13) 
(14) 
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From Equation (14), it can be seen that x is group sparse since Re(a) andlm(«) show the same 

sparse support, i.e., they are nonzero with the same indices, which would improve the recover 
performance if properly used [19]. Actually x has N groups corresponding to N complex range 

resolution cells with their real parts and imaginary parts as group members. By exploiting the group 
sparse a prior with ^ norm [20], Equations (11) and (12) can be transformed into: 

y=Ax+e (15) 
X = arg min ^ \\Ax - jI' + + (16) 

where x. is the /th element of x.With x estimated by Equation (16), the HRRP a can be achieved 
using Equation (14). N is the total number of the sparse groups, i.e., the total range resolution cells. 

3.2. Complex-Valued Fast Sequential Homotopy ( CV-FSH) Algorithm 

In this Subsection, we will extend the FSH approach to the complex domain. Assume that we have 
solved Equation (16) to get the estimation Xq as well as (a^ = Xq[1: N]-\- J-XqIN + I: 2N] ) with our 
current set of measurements for some given value of X . Now suppose that we get a new measurement 
given as si = u'a + n , where u' g C^""^ is a row vector and n denotes the complex noise in the new 
observation. Similarly, the real- valued model can be written as: 

"Re(^:) 
_Im(^;) 

which can also be regarded as adding two real-valued measurements sequentially. 

Herein, we consider adding the real part of s[ .The same way can be used to deal with the imaginary 
part. By defining w = Re(^') , b = [RQ(u'),-Jm(u')] and d = Rc(n) , the system of equations becomes: 



y 




A 




e 








x + 




w 




b 




d 



(18) 



Now we need to solve the following problem for update: 

x' = arg min -(||Ax - + \bx -wfj + ij^^xf + x^^ (19) 

Since the method used here is similar to the one in [14], we provide a description paralleling that of [14]: 

(1) Homotopy parameters establishment 

According to the homotopy technique, we need first to introduce continuous homotopy parameters 
linking the solved optimization program to the new one so that we can trace the solution path by varying 
the homotopy parameters carefully. From Equations (16) and (19), by introducing the parameter s , the 
new measurement is incorporated gradually as: 

X* = arg min - -yf^ + s \bx - ) + ^xf + xf^^ (20) 



Re(M'),-Im(M') 
Im(M'),Re(M') 





'Re(n')' 


x + 





(17) 
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While 6 increases from 0 to 1, the solved problem in Equation (16) changes to the new one in 
Equation (19). By resorting to the homotopy technique, we can trace the path of solutions, Le., from x 
to , as we slowly vary the parameter s . 

(2) Optimality Conditions Determination 

By differentiating the cost function in Equation (20) with respect to x , the optimality conditions for 
any solution x* to Equation (20) can be formulated as: 

Al^Ax^ -y) + sbl{bx^ -w) = -lD{xl) (21) 
\a^, (Ax* - j) + {bx - w)| < ^ |/)(x;, )| (22) 

where Fis the support of the solution x* , i.e., the nonzero or active element index set. is the 
complemental set. {H)^ denotes the new matrix/vector extracted from matrix/vector H with 

columns/elements indexed by F : 

'^o)/>/(^(o)'+(Wa^)'' r(0<A^ 

< (23) 

)\T{i)>N 

where |F| is the cardinality of the index set F , [-f denotes the transpose of a vector and F(/) denotes 
the /th element of F . 

Given the sparse support F , Equation (21) is the standard optimality conditions of nonlinear 
programming. Equation (22) holds since the sparsity exploitation term, the last in Equation (19), is more 
sensitive to the deviation of inactive elements than the representation error term shown in the bracket, 
i.e., the inactive elements contribute mainly for the sparsity enhancement and little for the representation 
error reduction. 

(3) Solution Path Tracing 

Assume we are at one of these critical values of s = s^ and the solution x^ has been achieved. 
Adding s with an infinitesimal amount to 8^' , we can evaluate the corresponding solution x^ with 
dx = x^ -x^ by using the optimality conditions Equation (21) and first-order Taylor approximation: 



D(x*) = [Ji,J2,...,^|r|f ' = 



dx = 



-(s; -8j(A,^A, +D\{x,\) + E;blb^ybl{bx, -w\ on F ^^^^ 
1 0, otherwise 



where: 



D'{{x^\) = diag{id^,d^,...,d^A), d. = 



The derivations of Equations (24) and (25) are given in the Appendix A. Through reformulating 
dx = Qj^dx , we can derive the update direction dx and the positive step size 9^ as follows using the 
matrix inversion lemma. The details are presented in the Appendix B: 

^ \-U-X(bx,-w% onF 

n u . (26) 

0, otherwise 
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(27) 



l + (s/ -8jw 

where U = A^A^ +DX(Xj^\) + s,^b^b^ , u = b^U'^b^ . 

The solution moves away from along the direction dx as s increases, until the current support T 
must be modified. That means certain existing element of indexed within F will shrink to zero or certain 
one indexed within will be active. For the last case, the corresponding inequality in Equation (22) will 
become an equation. Furthermore, the modification is made according to the required smallest step size. 

For the active element shrinking to zero, i.e., turning into the inactive one, the smallest step size is: 



9 = min 

Yer,Y<A/' 



dx(y) 



+ 



dx(y + N) 



(28) 



where y is the group index. Since Xj^_^^(y) = Xj^(y) + Qdx(y) do certainly not shrink toward zero with 



negative 



dx(y) 



and positive step size, we use (•)^to denote that the minimum is only taken over the 



group indices with both parts shrinking toward zero simultaneously. Assuming the group y needs the 
smallest step size. 

For the inactive element turning into the active one, the required step size is what makes any one of 
the inequalities in Equation (22) become an equation, which can be written as: 



A^, (Ax,' - J) + s.'blc (bx^ -w)<X Diix^ \, ) 



(29) 



where the equation holds for some j e F^ , which means the inactive element j is turning into the 
active one. 

Substituting x^' = + ^k^^ ii^to Equation (29) and setting: 

= A" (Axj^ - j) + Sj^b^ (bXj^ - w) 



(30) 

dj^ =(A'A-\-Sj^yb)dx+y(bXj^-w) (31) 

Then Equation (29) can be rewritten as: 

(p,).+Q,(d,).^±W({x;).) = ±W{(x,).) (32) 

where the second equation holds since (dx)j = 0 . The first approximation is deduced in Appendix C. It 
should be noted that Equations (30) and (31) are consistent with Equation (35) in [14]. Thus, the smallest 9^ 
satisfying Equation (32) is given by: 

9^= min i9(y)f +[9(y + A^)f 



(33) 



where: 



^(y) = min 



A (Y) = X, (Y)/ ^[x,(y)f+[x,(y + N)f 



(34) 
(35) 
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9(y) is the smaller step size to achieve Equation (32) for the real part of group y , with %{y + N) for 
imaginary part. Suppose that the group needs the smallest step size. 
Thus the step size corresponding to the critical point is: 

e=min(e^e-) 06) 

If 6 = 9^, the corresponding inactive element turns into the active one. Otherwise, the corresponding 
active element turns into the inactive one. 

Then, according to Equation (27), the new critical point is: 

9 



and the updated solution is: 

Based on the above idea, CV-FSH can be described by Algorithm 2. 



(37) 



(38) 



Algorithm 2 Dynamic update with sequential complex measurements 

Step 1: Initialize as solution to Equation (16) with support T . 

Step 2: Adding the real part of the measurement. 

Set ^ = 0, So=0, w = Re(5;) and ft = [Re(ii'),-Im(ii')] . 

Step 3: Iteration: 

(1) Compute the update direction dx with Equation (26). 

(2) Compute the step size 9 with Equations (28), (33) and (36). 

(3) Update critical point and the solution with Equations (37) and (38), respectively 



(4) If e... 
then e = 



>1 

1-8,. 



= +95x, s^+i =1 



l + (l-sji/ 

break (quit the loop); 
end if 

(5) If9 = 9- 
r^r\{y",y"+A^} 

else 

r^rU{y",y"+A^} 

end if 

(6) k^k + \ until stopping criterion is satisfied. 
Step 4: Adding the imaginary part of the measurement. 

A 

Set A= ^ , , v= _ " . , w 





y 


, y = 


Re(<)_ 



[Re(«'),-Im(M')] 
Sq = 0 , and execute the Iteration (Step 3) again. 
Output: The updated solution f . 



= \m{s[ ) ,b = [lm(M'), Re(«')] = x,,, , k = 0 , 
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4. Stopping Rules 

When the number of pulses is not sufficient, the results of the sparse recovery are inaccurate and will 
change notably even when only one new measurement is added. Otherwise, the results are accurate and 
become stable with new added measurements. Therefore, we can define the stopping rules mainly from 
the one-step agreement of two sequential results. 

Obviously, one of the necessary conditions to stop the processing is that the radial velocities are 
equivalent for the two sequential processing results, i.e., =y? ^ which will guarantee that the 

motion compensation is well performed. This is because the estimated mean velocity remains nearly 
unchanged if only one more PRI is added to CPI, which is comparatively much longer. Another stopping 
criterion is the sufficient similarity of the HRRPs obtained by the two sequential processing. Assume 
that for a certain instance the minimum number of required pulses is M^^^ which just satisfies the RIP 

condition, then we can expect that, when M < M^^^ , the similarity between two adjacent reconstructed 

HRRPs will be relatively low because of the large error in CS reconstruction, whereas the two adjacent 
reconstructed HRRPs will be matched very well when M > M^^^ . Therefore, the value of M^^^ can be 

estimated by measuring the normalized correlation coefficient between the Mth and (M + l)th HRRPs, 
which can be given as [21]: 

Cm^i= u. II ||. II (39) 

IrM II2 |nM+l||2 

If the correlation coefficient is smaller than a certain threshold c^^ , the value of M must be increased 

and one more pulse whose carrier has not been used ever should be transmitted accordingly. On the 
contrary, a correlation coefficient larger than c^^ implies that the HRRP has already been correctly 

generated and the dynamic processing can be stopped. 

Finally, it can be concluded that for the proposed dynamic algorithm, if the number of pulses satisfies 
the following three conditions simultaneously, the HRRP and estimated velocity of the target can be well 
obtained with high probability: (1) The number of pulses is among the presetting scope; (2) The 
estimated velocity is equivalent with the previous one; (3) The normalized correlation coefficient 
between the current HRRP and the previous one is sufficiently high to approach one. 

5. Experiments 

In this section, the feasibility and performance of the proposed dynamic HRRP generation algorithm 
are tested on both simulated and chamber measured data. Our simulations are performed in the 
MATLAB7 environment using a Pentium (R) 4 CPU 3.00 GHz processor with 1 GB of memory. 

5.1. Effectiveness of Complex-Valued FSH ( CV-FSH) 

In this subsection, simulations with an ideal scattering center model are carried out to validate the 
effectiveness of the complex- valued FSH algorithm. The parameters of the transmitted RSFR signal are 
set as follows: the carrier frequency ranges from 9.5 GHz to 10.5 GHz, Tr =100 |lis and r= 0.1 |lis. It can 
be calculated that = 1 ns and N = 100. Assume that the numbers of targets' scattering centers K are 
set to be 20, 30, respectively and they are located within length 15 m in the LOS. The real parts and 
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imaginary parts of their amplitudes are both randomly distributed A/'CO,!^) . The distance from the target 
to the radar is 6 km and the target is assumed to be stationary. The minimum and maximum numbers of 
transmitted pulses are 30 and 90, respectively. When the target has been detected, an RSFR signal is 
transmitted to obtain the HRRPs of the target under white Gaussian noise, the real parts and imaginary 
parts of which are both randomly distributed MiO.OSdl^) . For comparison, the complex- valued FSH 

(CV-FSH) algorithm as well as real-valued FSH (RV-FSH) algorithm and the well-known Convex (CVX) 
toolbox [22] are applied to 200 Monte Carlo trials. For RV-FSH, complex-valued sparse recovery 
Equation (11) is transformed into the real- valued one Equation (15). However, RV-FSH exploits the 
conventional sparse a prior only, compared to CV-FSH imposing group sparse a prior. Thus, RV-FSH is 
a direct application of the approach proposed by [14]. 

The average normalized correlation coefficients (Ancc) between the actual scene and the recovered 
scenes are calculated for each step and the results are shown in Figure 3a, from which it can be seen that 
with all three algorithms the Ancc rises quickly with the increase of the number of pulses and then tends 
to a constant approaching one when the number of pulses is large enough. Furthermore, it can also be 
seen that with the increase of the number of pulses and the decrease of the number of scattering centers, 
e.g., M > 3K , the RIP condition is gradually satisfied, which will lead to accurate parameter estimation 
with the recovery algorithm. For the precise solution algorithm of convex optimization problem 
Equations (12), (16) and (19) in the CVX toolbox, it performs slightly better than our proposed CV-FSH 
algorithm, which is a fast pursuit-like method to the update problem (19) and approaches the CVX 
solutions quickly as the number of pulses increases. RV-FSH performs worst since it imposes a weak 
prior or loses some information from the information theory perspective. 

Figure 3b presents the average CPU time required to complete each algorithm. The curves show that 
the cost time of the CVX toolbox increases for larger scale of the sensing matrix. However, for the 
RV-FSH and the CV-FSH algorithms, the cost times almost remain unchanged and even decrease when 
the results are stable. This is because when the results change slightly toward the true solution, it can 
provide better and better initialization for faster converging in the updating procedure. In Figure 3b, it 
can be seen that the CPU times required for RV-FSH and the CV-FSH algorithms begin to reduce notably 
when the number of pulses is about three times of the number of scattering centers. This phenomenon occurs 
since the support remains almost unchanged when the number of pulses overpasses certain number and so 
only few iterations are required to update the solution. CV-FSH achieves the stable point more quickly 
where the cost time begins to be low and almost unchanged owing to its superior recovery performance. 

Therefore, considering both computational complexity and the recovery performance, we observe 
that the CV-FSH algorithm is more suitable for our problem. It is a great potential to update the solution 
real-time or quasi-real-time with sequential measurements by resorting to CV-FSH. 

The CV-FSH algorithm is further analyzed concerning its robustness of HRRP generation. Here we 
mainly consider how the algorithm performs when adding different white Gaussian noises. By setting 
the variance of noise as 0.01, 0.1, 1 (corresponding to SNR 20 dB, 10 dB and 0 dB), respectively, 100 
Monte Carlo simulations are performed for each case and the Ancc between the actual scene and the 
recovered scenes are calculated. As shown in Figure 4, as the SNR decreases, the Ancc converges more 
slowly, i.e., the required number of measurements for exact recovery increases. When the SNR is too low 
{i.e., SNR = 0 dB), the algorithm is almost ineffective. This phenomenon is similar to the CVX toolbox [22]. 
Actually, the descending performance of sparse-based algorithms under lower SNR is still an open 



Sensors 2014, 14 



8296 



problem. In the actual radar imaging applications, the required SNR of echo data for imaging is often 
higher than that for detecting and tracking. Therefore, from the engineering point of view, the stability of 
the algorithm with SNR higher than 10 dB can achieve most application requirements. 



Figure 3. HRRP generation performance of three algorithms with various number of 
scattering centers, (a) Average correlation coefficients between the actual scene and the 
recovered scenes versus number of pulses; and (b) Average CPU time required to complete 
each algorithm versus number of pulses. 
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Figure 4. HRRP generation performance with various numbers of pulses under different SNRs. 
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5.2. Sequential HRRP Generation of Stationary Target 



Before testing the feasibility of Algorithm 1, we first make some experiments to validate the 
effectiveness of the SCS algorithm for sequential HRRP generation of the stationary target in Figure 5a. 
The radar for anechoic chamber measuring works in the sweep mode, whose carrier frequencies change 
from 9 GHz to 11 GHz with an interval of 20 MHz and the center frequency is =10 GHz. The 

elevation angle of the target is 0°, and the azimuth angle varies from 0° to 90° (nose direction is 0°) with an 
interval of 0.2°. By performing 128-point sparse recovery {i.e., the number of range cells is set to be 128) on 



Sensors 2014, 14 



8297 



the 100 samples of the wideband echo signal for each azimuth cell, we can obtain all 450 HRRPs of the 
target, which are shown in Figure 5b. It should be noted that in the following analysis the HRRP 
obtained from the whole raw data are called the original HRRP. 



Figure 5. (a) Model of the anechoic chamber measured target; and (b) Normalized HRRPs 
for all azimuths. 
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As the RSFR is still under development, for convenience, herein we randomly draw out M (from 16 to 64) 
measurements from 100 samples of the LSFR to simulate the echo signal for RSFR. For stationary 
targets, this is equivalent to the samplings of output signal in RSFR which transmits the M pulses at 
corresponding frequencies. Since the selection of the samples intervening in the SCS approach is done 
randomly, the repeatability of the proposed method has to be checked. In order to do so, 200 Monte 
Carlo simulations are performed at two typical azimuths (5° and 45°) to determine the proper number of 
measurements. In a general way, choosing a larger number of experiments is more reasonable. However, 
an appropriately small number is preferred for the low computation load. What is more, hundreds or 
even tens of independent trials are sufficient in the sense that the outcome change slightly when more 
independent trials are added, which is observed in our simulations. Nevertheless, strict analysis 
especially mathematical justification about the appropriate number of experiments is still valuable. This 
open problem deserves further research. 

Firstly, the original HRRPs for the two azimuths are given in Figure 6a,b, which obviously show that 
the number of strong scattering centers for azimuth 45° is smaller than that for azimuth 5°. Due to its 
much larger computation load, CVX is not appropriate at all to implement the fast sequential update, the 
core of our dynamic compressive sensing strategy. Thus, the CV-FSH approach is applied to the selected 
measurements and compared with RV-FSH. Figure 6c and d present the Ancc between two adjacent 
HRRPs (Ancci) for various numbers of measurements. For comparison, the Ancc between the recovered 
HRRP and the original HRRP (Ancc2) are also given. For CV-FSH, it can be seen that when the HRRP 
tends to be stable at about 44 measurements for azimuth 5° {i.e., Ancc2 is larger than 0.95), the Ancci is 
always larger than 0.99. Therefore, we can regard 44 as the most proper number of measurements. In the 
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case of azimuth 45° with fewer scattering centers, similar phenomenon can be found and the most proper 
number of measurements is about 21. 



Figure 6. Results of the SCS algorithm for HRRP generation of the stationary target, 
(a) Original HRRP for azimuth 5°; (b) Original HRRP for azimuth 45°; (c) Ancci and Ancc2 
for azimuth 5°; and (d) Ancci and Ancc2 for azimuth 45°. 
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For RV-FSH, the Anccs are lower than those of CV-FSH. The most proper number of measurements 
is 50 for azimuth 5° and 24 for azimuth 45°, both of which are larger than their counterpart of CV-FSH. 
These results mean that our CV-FSH outperform RV-FSH on reducing the required number of pulses for 
HRRP synthesis. This experiment also tells that as long as the number of measurements is sufficiently large, 
precise HRRPs with slightly decreasing performance can be obtained with much fewer measurements than 
traditional sense. It is reasonable because the CS-based algorithms perform stably as far as the RIP condition 
is satisfied, but much worse if the condition is not well satisfied. 

In order to further validate the effectiveness of the stopping rules, the proper numbers of 
measurements for generating high quality HRRPs at all azimuths (0° to 90° with an interval of 1°) are 
calculated. The stopping rules are that the number of measurements is among the presetting scope (from 
16 to 80) and the Ancci is larger than 0.99. The original HRRPs are shown in Figure 7a and the proper 
numbers are presented in Figure 7b, from which we can see that with the reduction of the number of 
strong scattering centers projected into the radar LOS, the required number of measurements also 
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decreases correspondingly. Therefore, we can conclude that the number of measurements should be 
nearly proportional with the number of scattering centers, which is also consistent with the RIP 
condition. It can be seen again that the required numbers of pulses when resorting to CV-FSH are 
consistently smaller than that using RV-FSH. 

Figure 7. Proper numbers of measurements for HRRP generation, (a) Central supports of 
normalized HRRPs for all azimuths; and (b) Proper numbers of measurements for all azimuths. 



50 



55 

</> 

1 60 
O 

0) 

P 65 

03 

70 
75 
80 



20 40 60 80 
Azimuth Cells 

(a) 

5.3. Sequential HRRP Generation of Moving Target 
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Finally, the performance of velocity measurement and HRRP generation with Algorithm 1 is verified. 
As the target used for chamber measuring is with no translational motion, we will simulate the echo data 
by assuming that the target is moving toward the radar with a certain radial velocity v and the selected 
measurements s^{m\ are multiplied by the synthesized phase terms corresponding to the radial velocity 

exp(j27rF^2vmr / c) , i.e.: 



s[ [m] = s^ [m] • exp( j27rF^ 2vmr / c) 



(40) 



We still perform the experiments for two typical azimuths (5° and 5°). In order to make the scenario 
more practical, we assume that the actual radial velocities of the target at azimuth 5° and 45° are 
respectively 25 m/s and 75 m/s, which is smaller than the maximum unambiguous velocity 
c I (2/^r ) = 150 m/s when the pulse repeat time is 100 us. For the target with higher radial velocity, e.g., 

200 m/s, we usually measure its ambiguous velocity as 50 m/s, which belongs to [0 m/s, 150 m/s] , rather 
than its true value. Therefore, here we select two typical values {i.e., 25 m/s and 75 m/s ) for the target 
velocities. Assume that the accuracy of tracking system is 5 m/s and measured velocities are 25 m/s 
and 75 m/s , respectively, then we can regard the velocity domains of interest as [20 m/s, 30 m/s] and 
[70 m/s, 80 m/s] . The velocity domain is divided into L = 10 points and then by applying Algorithm 1 
to the new data set, we can obtain the velocity estimation as well as the HRRPs. 

In order to see the statistic behavior, the experiment is also performed based on 200 Monte Carlo 
simulations. As shown in Figure 8a, when the number of measurements is larger than 39 for CV-FSH 
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approach or 40 for RV-FSH approach, the velocity estimation of the target at azimuth 5° will be always 
precise, and for azimuth 45° the minimum number of measurements is about 23 for CV-FSH approach 
or 26 for RV-FSH approach. Herein, apart from the precise velocity, the Ancc between two adjacent 
HRRPs should also be large enough, which will finally determine the optimal number of pulses required 
in each CPI. As shown in Figure 8b, the values of Ancc between two adjacent HRRPs are always larger 
than 0.99 when the number of measurements exceed 47 for CV-FSH approach or 49 for RV-FSH 
approach at azimuth 5° and 25 for CV-FSH approach or 28 for RV-FSH approach at azimuth 45°, which 
are close to the results for the stationary situation in Figure 6. It can be seen again that our CV-FSH 
algorithm outperforms the RV-FSH algorithm. 

Figure 8. Results of the SCS algorithm for velocity measurement and HRRP generation of 
the moving target, (a) Velocity estimations for various numbers of measurements; and (b) 
Ancc between two adjacent HRRPs. 
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6. Conclusions 



(1) A dynamic compressed sensing strategy for HRRP generation using RSFR is presented in this 
paper. In this strategy, a complex- valued fast sequential homotopy (CV-FSH) algorithm is proposed to 
implement the HRRP update quickly with sequentially transmitted and received pulses. By transforming 
the complex- valued SCS signal model to the group sparse real-valued one a priori and adding the real 
part and imaginary part one by one, CV-FSH performs as an efficient recursive sparse recovery 
procedure. The stopping rules for this dynamic strategy are based on the one-step agreement of velocity 
and HRRP. If the sequential estimations of target redial velocity are equivalent and the normalized 
correlation coefficients between two adjacent HRRPs are sufficiently large, the minimum number of 
pulses required in each CPI can be determined. The results using simulated and real data show that the 
proposed algorithm is more suitable for practical application and especially for the uncooperative targets. 

(2) Algorithm for fast sequential update is the core of our dynamic compressive sensing strategy. 
Low computation complexity is the main requirement for this algorithm. It has been demonstrated that 
our CV-FSH algorithm makes a better tradeoff between reconstruction accuracy and computation 
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complexity, thus it is better for fast sequential update. In the aspect of reconstruction accuracy, 
simulation results show that CV-FSH is slightly inferior to the algorithm in CVX tool box for precise 
solution. However, computation load of the CV-FSH is almost 100 times lower than that of the 
algorithm in CVX tool box under our simulation settings. Moreover, when more and more measurements 
are available, the computation complexity remains unchanged and even decrease for the CV-FSH 
algorithm since the results become more and more stable and so fewer iterations are required, while the 
computation complexity of the algorithm in CVX tool box increases for the larger scale of the sensing 
matrix. These mean that the algorithm in CVX tool box is not appropriate at all to implement fast sequential 
update. Compared with RV-FSH algorithm (a direct application of the approach proposed by [14]), 
CV-FSH algorithm performs better both on reconstruction accuracy and computation complexity. The 
results using real data show that the required number of pulses when resorting to CV-FSH is reduced 
consistently by about 10%. 

(3) Due to the efficiency and universality, the CV-FSH algorithm can also be used in other CS-based 
radar applications such as moving target detection [23,24] and cross-range compression in 
two-dimensional imaging [25,26]. Furthermore, extending the dynamic HRRP generation to the 
polarimetric case is somewhat straightforward but significant. 
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In this appendix, we will give the deviation dx of the solutions arising from homotopy parameter 
increase using the optimality conditions. The optimality conditions are given: 



Subtract the left terms of Equation (41) from the left of Equation (42) and make a first-order Taylor 
approximation to achieve: 



AlA{x^ -x^) + E^blb{x^ + -^k)blbx^+{z^ -z^)blw+W\{x^\){x^ - xj « 0 (43) 



and so the deviation dx can be obtained by simplifying the equations in Equation (43). 
Appendix B 

In this appendix, the update direction of the solutions and the step size are derived in details. These 
are separated in the deviation dx of the solutions, which is simplified as: 



Appendix A 



Al {Ax, -y) + E,b^ (bx, - w) + W((x, )^) = 0 



(41) 



A^ (Ax,' - j) + s,'b^ (bx^ - w) + XD((x/ = 0 



(42) 



^x = -^^{AlA^ +D\{x,\) + E^blb^y'bl(bx,-w) 

= -As( + DX(x, )r ) + s,b^b^ + Asb^b^ y'bl (bx, - w) 



= -As[U-' -U-'bl{— + uy\U-']bl{bx,-w) 



As 



= -Ae[U-X (bx, - w) -U-'bl ibx, - w)] 



(44) 



= -A8[l - ( — + u) ^u\U ^bl (bx, - w) 



Thus the step size is: 



Q,=Ae[l-(— + uy'u] = 



1 + Asw 



(45) 



and the direction is given as: 




(46) 



Appendix C 

In Equation (29), the inactive element j is turning into the active one, which means: 
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Aj {Ax^ -y) + e/Z»J (bx^' -w) = X 

substituting x^.' = aTj, + ^(^Sr into Equation (47): 

A] (Ax,' -y) + e,'b] (bx,' - w) 
= A](Ax,-y) + e,b](bx,-w) + dMjA + e,b]b)dx + (e; -e,)b] (bx,+d,bdx -w) 

J J J (48) 

= A] {Ax, -y) + E,b] (bx, - w) + 0, (Aj A + E,b]b)dx + (9, + 9, A8,w)ftJ (bx, + Q.bdx - w) 
« Aj (Ax, - J) + 8,6J (6x, - w) + 9, (Aj A + E,b]b)dx + 9,6j {bx, - w) 

where As, = (s,' -s,) . Since As, <^ land 9, <As, , the approximation is reasonable. Thus, the 
approximation in Equation (32) holds. 
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